Local regression methods and systems for image processing systems

ABSTRACT

This disclosure provides methods, apparatus and systems for performing image processing regression for approximating multidimensional color transformation. According to an exemplary method, a shaping matrix is selected to minimize a cost function associated with a local linear regression representation of the color transformation. In addition, an alternating least squares algorithm is utilized to jointly optimize regression and shaping parameters.

CROSS REFERENCE TO RELATED PATENTS AND APPLICATIONS

U.S. patent application Ser. No. 12/406,303, filed Mar. 18, 2009, entitled “LOCAL REGRESSION METHODS AND SYSTEMS FOR IMAGE PROCESSING SYSTEMS,” by Monga et al. is incorporated herein in its entirety.

BACKGROUND

Local linear regression is used in a variety of data fitting applications. Particular applications within the realm of color imaging include printer and scanner characterization. A typical regression problem involves first gathering a training set of input data points from an input space and corresponding output data points from an output space. For the color characterization application, both input and output spaces are multi-dimensional color spaces. The goal of the regression algorithm is then to derive mappings from every point in the input space to the output space while minimizing error over the training set. An additional consideration is to ensure that the regression does not overfit the data in the sense that it is robust enough to filter out noise in the training data. Local regression algorithms are often used in situations where a single global fit may be inadequate to approximate complex non-linear transforms, as is typical in printer characterization. Instead, local transforms are derived where the regression parameters vary as a function of the input data point. Locality in regression is achieved by using a weighting in the error minimization function which varies (typically decays) as a function of the distance from the regression data point. Choice of these weight functions is typically intuitively inspired, and not optimized for the training set. This sometimes results in large regression errors especially with sparse training data. A key fundamental question hence remains on how to best use a certain local neighborhood of data points in regression problems and how to quickly compute optimally local transforms using certain local neighborhoods of data points in regression problems.

Incorporation By Reference

The following references are totally incorporated herein by reference.

R. Bala, “Device Characterization,” Digital Color Imaging Handbook, Chapter 5, CRC Press, 2003.

C. G. Atkeson, A. W. Moore and S. Schaal, “Locally weighted learning,” Artificial Intelligence Review, 1997.

E. Chong, S. Zak, “An Introduction to Optimization,” 2^(nd) Ed, Wiley, 2001.

BRIEF DESCRIPTION

In one aspect of this disclosure, described is a regression method for approximating a multidimensional color transformation comprising (a) receiving a set Γ of training samples (x_(i), y_(i)), 1≦i≦T, where x_(i) represents input color data to the multidimensional color transformation, and y_(i) represents corresponding output color data from the multidimensional color transformation; (b) receiving an input color x; (c) selecting a regression function f(x) parameterized by a regression matrix A_(x) that approximates the multidimensional color transformation at the input color x; (d) generating a cost function C representing a localized error produced by the regression function f(x) on the training set Γ, where the localized error is a function of both the parameters of f(x) and a shaping function that defines the shape and orientation of a neighborhood of training data localized around the input color x; (e) deriving the parameters of the regression function f(x) and shaping function to jointly minimize the cost function C by using an iterative alternating least squares algorithm to determine the elements of the regression matrix A_(x) and the shaping function; and (f) generating an output color y by calculating f(x) using the derived parameters of step (e).

In another aspect of this disclosure, described is a computer program product, that when executed by a computer, causes the computer to execute a regression method for approximating a multidimensional color transformation, the method comprising (a) receiving a set ″ of training samples (x_(i), y_(i)), 1≦i≦T, where x_(i) represents input color data to the multidimensional color transformation, and y_(i) represents corresponding output color data from the multidimensional color transformation; (b) receiving an input color x; (c) selecting a regression function f(x) parameterized by a regression matrix A_(x) that approximates the multidimensional color transformation at the input color x; (d) generating a cost function C representing a localized error produced by the regression function f(x) on the training set Γ, where the localized error is a function of both the parameters of f(x) and a shaping function that defines the shape and orientation of a neighborhood of training data localized around the input color x; (e) deriving the parameters of the regression function f(x) and shaping function to jointly minimize the cost function C by using an iterative alternating least squares algorithm to determine the elements of the regression matrix A_(x) and the shaping function; and (f) generating an output color y by calculating f(x) using the derived parameters of step (e).

In still another aspect of this disclosure, described is an image processing method for rendering an image on an image output device comprising receiving a device independent color space representation of the image; accessing an inverse characterization transform associated with the image output device to generate a device dependent color space representation of the image, the inverse characterization transform representing the inverse of a multidimensional color transformation associating a plurality of device dependent color space values with a plurality of respective device independent color space values, the multidimensional color transformation generated by performing a method comprising (a) receiving a set Γ of training samples (x_(i), y_(i)), 1≦i≦T, where x_(i) represents input color data to the multidimensional color transformation, and y_(i) represents corresponding output color data from the multidimensional color transformation; (b) receiving an input color x; (c) selecting a regression function f(x) parameterized by a regression matrix A_(x) that approximates the multidimensional color transformation at the input color x; (d) generating a cost function C representing a localized error produced by the regression function f(x) on the training set Γ, where the localized error is a function of both the parameters of f(x) and a shaping function that defines the shape and orientation of a neighborhood of training data localized around the input color x; (e) deriving the parameters of the regression function f(x) and shaping function to jointly minimize the cost function C by using an iterative alternating least squares algorithm to determine the elements of the regression matrix A_(x) and the shaping function; and (f) generating an output color y by calculating f(x) using the derived parameters of step (e).

In yet another aspect of this disclosure, described is a computer program product that when executed by a computer, causes the computer to perform a color transformation for rendering an image on an image output device, the method of deriving the color transformation comprising (a) receiving a set Γ of training samples (x_(i), y_(i)), 1≦i≦T, where x_(i) represents input color data to the multidimensional color transformation, and y_(i) represents corresponding output color data from the multidimensional color transformation; (b) receiving an input color x; (c) selecting a regression function f(x) parameterized by a regression matrix A_(x) that approximates the multidimensional color transformation at the input color x; (d) generating a cost function C representing a localized error produced by the regression function f(x) on the training set Γ, where the localized error is a function of both the parameters of f(x) and a shaping function that defines the shape and orientation of a neighborhood of training data localized around the input color x; (e) deriving the parameters of the regression function f(x) and shaping function to jointly minimize the cost function C by using an iterative alternating least squares algorithm to determine the elements of the regression matrix A_(x) and the shaping function; and (f) generating an output color y by calculating f(x) using the derived parameters of step (e).

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 is a flow chart of an exemplary image processing regression method according to this disclosure.

FIG. 2 is a plot of weight functions w(x, x₀) for x₀=(0,0), where the weight function decays as a function of distance from x₀ to generate a neighborhood around x₀.

FIG. 3 illustrates one exemplary function mapping from 2-D space (R²) to 1-D space (R) to be approximated according to an exemplary regression method according to this disclosure.

FIG. 4 illustrates the contours of the function illustrated in FIG. 3 with exemplary training data overlaid.

FIG. 5 illustrates local regression using neighborhood shaping according to an exemplary embodiment of this disclosure.

FIG. 6 is a block diagram of an image processing system using an exemplary regression method according to this disclosure.

FIG. 7 is a flow chart of another exemplary image processing method according to this disclosure.

DETAILED DESCRIPTION

This disclosure provides methods and systems for local regression in deriving color transformations by introducing the notion of “shaping” in the localizing weight function. The disclosed exemplary embodiments include: 1.) a parameterization of the weight function typically used in local regression problems via a shaping matrix, and 2.) a method to obtain the “optimal” shaping matrix by explicitly introducing the weight function parameters in the regression error measure. Demonstrated experimentally are that significant gains can be made by optimizing “the shaping matrix” in local regression problems. Many color imaging applications including printer and scanner characterization can benefit from the disclosed methods, apparatus and systems. The disclosed exemplary embodiments are particularly advantageous for color devices that employ a large number of color channels, thus inducing a large dimensionality in the characterization data.

In addition, this disclosure provides practical methods and systems of computing the optimal shape parameters in local regression via a relatively fast algorithm. This makes the disclosed regression techniques applicable to large and high-dimensional data sets. The fundamental basis of the algorithm is to identify a weighting/shaping function for which the regression cost function is separably convex. That is, though not jointly convex in the regression and shape parameters, the cost function is convex with respect to regression parameters when the shape parameters are fixed, and conversely, the cost function is convex with respect to shape parameters when the regression parameters are fixed. This allows the development of a fast algorithm based on an alternating least squares (ALS) technique.

Benefits are apparent over standard local regression and/or the use of naïve gradient methods to optimize the parameters jointly.

Regression is a common technique for estimating a functional relationship between input and output data, and is used frequently to derive color device characterization transformations. The latter typically establish a functional relationship between a device dependent color space and a device independent color space. For printers, examples of device dependent color spaces include CMY, CMYK or CMYKOV, where the symbols stand for Cyan, Magenta, Yellow, Black, Orange, Violet, respectively. For display devices, the prevalent device dependent color space is RGB (or Red, Green, Blue). A common example of a device independent color space is CIELAB. There are two types of color characterization transforms—a forward and an inverse. For output devices, the forward transform maps a device dependent color space to a device independent color space, and conversely, the inverse transform maps a device independent color space to a device dependent color space. In certain applications, the forward transform for one device is concatenated with the inverse transform for another device to produce a “device-to-device” color transformation. The regression techniques and exemplary embodiments described herein can be applied to forward, inverse, or device-to-device characterization transforms.

Linear regression is a specific case where the functional relationship between the input and output spaces is approximated by a linear transform. When the input and output data belong to multidimensional vector spaces, the linear transform is a matrix. Specifically, consider the problem where y in R^(m) is to be estimated as a function of an input variable x in R^(n). Thus we have y≈f(x). Let Γ={(x_(i), y_(i)), i=1, 2, . . . , T} denote the set of training data over which this response in known. The linear approximation is given by: y=ƒ(x)=A·x, xεR ^(n) , yεR ^(m) , AεR ^(m×n).  (1)

The “best” regression parameter A is determined by minimizing the regression cost function that describes an aggregate error between y_(i) and A_(xi) for the training set.

A variant of this is local linear regression as described in R. Bala, “Device Characterization,” Digital Color Imaging Handbook, Chapter 5. CRC Press, 2003, wherein the matrix A varies as a function of location x in input space. Thus, we have: y=ƒ(x)=A _(X) ·x, xεR ^(n) , yεR ^(m) , AεR ^(m×n).  (2)

For each input data point x, the “best” regression parameter A_(x) is determined by minimizing the regression cost function:

$\begin{matrix} {{C\left( A_{x} \right)} = {\frac{1}{T}{\sum\limits_{i = 1}^{T}{{{y_{i} - {A_{x} \cdot x_{i}}}}^{2} \cdot {{w\left( {x,x_{i}} \right)}.}}}}} & (3) \end{matrix}$

In the above cost function, note that it is the presence of weight function w(x, x_(i)) that introduces locality in regression.

The most general requirement is for w(x, x_(i)) to decay as a function of the distance d(x, x_(i))=∥x−x_(i)∥. A popular instantiation is: w(x,x _(i))=e ^(−α(∥x−x) ^(i) ^(∥) ² ⁾  (4)

The above weight function is plotted in FIG. 2 for a 2-D input variable x. Clearly, there is decay as a function of distance from input point x which means that in the cost function of Equation (3), more weight is attached to regression sets for which d(x, x_(i)) is small. Qualitatively, the subset of {x_(i)}'s for which w(x, x_(i)) is greater than a threshold (so as to have a significant impact on the cost function) constitutes a local neighborhood of x.

Limitations of existing local linear regression is now discussed.

The use of such locality inducing weight functions is well-known to significantly help with regression accuracy over using a single “global” regression. The same has also been successfully applied to derive both forward and inverse printer color transforms in R. Bala, “Device Characterization”, Digital Color Imaging Handbook, Chapter 5, CRC Press, 2003.

The notion of locality as in Equation (4) is meaningful from the viewpoint of neighborhood size, i.e. a certain α may be chosen to control the spread of w(x, x_(i)) around x. That said, an important consideration that was previously ignored is shaping of w(x, x_(i)).

To appreciate this, consider, the 2-D function plotted in FIG. 3. FIG. 4 shows contours of this function with training data overlaid. It is desired to approximate this 2-D function with a locally linear regression function f( ) and to compute the regression output for the input point labeled “x”.

It may be seen from FIG. 4 that many different regression neighborhoods may be defined in the vicinity of the input point. A few such circular neighborhoods are shown in FIG. 4. Note also that circles centered at a point x₀ are contours of the distance function d(x, x₀)=∥x−x₀∥.

From a visual inspection of the plot in FIG. 3 it is clear that the original 2-D function is simply linear in the vicinity of the input point x, and can hence be perfectly approximated. But for that to happen, the regression must pick or attach greater weight to training data in the “linear region” vs. training data points in the non-linear region. Note further that for the three different circular neighborhoods shown in FIG. 4 the training data points identified as “x”s are closer to the input point, and will in fact get more weight, regardless of the size (radius) of the circular neighborhoods.

Detailed next is how this problem may be averted by introducing the notion of “neighborhood shaping.”

In this disclosure, utilized is the use of neighborhood shaping to address the problem discussed above. That is, the crucial observation is made that for a fixed neighborhood size, shaping plays an important role in regression accuracy.

Described now is how to achieve a desired neighborhood shaping. It is observed that for finite-dimensional input spaces, the distance function d(x, x_(i))=∥x−x_(i)∥ can alternatively be written as ∥x−x _(i)∥²=(x−x _(i))^(T)(x−x _(i))

As shown in FIG. 4, contours of such a distance function result in hyper-spheres in R^(n) (special case circle in R²).

Proposed first is a generalization of this distance to: ∥x−x _(i)∥_(Λ)=(x−x _(i))^(T)Λ(x−x _(i)) where Λ is a positive definite matrix, which is a requirement to ensure non-negativity of the distance metric for all x, x₀.

It is clear now that the contours of this new distance can be generalized to be elliptical. A diagonal Λ with positive unequal diagonal entries results in a hyper-ellipse with different ellipse radii in different dimensions, while non-diagonal choices of Λ allow the control of orientation.

Notably, the local linear transform A_(x) and the resulting output estimates may vary considerably with different choices of Λ. One possible strategy to optimize Λ is to make it proportional to the sample covariance matrix of the training data.

This disclosure and the exemplary embodiments described herein provide methods to find the best “shape” of the weight function or equivalently Λ for a fixed size/volume. To formally distinguish shape from size, Λ is re-written as follows: Λ=λS^(T)S; Λ, SεR^(n×n), where S denotes a “shape matrix” with determinant 1, and λ is a non-negative scalar relating to the size of the local neighborhood.

Given this separation of size and shape the shaping matrix may be solved for optimally by minimizing the regression error:

$\begin{matrix} {{{{C\left( {A_{x},\Lambda} \right)} = {\frac{1}{T}{\sum\limits_{i = 1}^{T}\left. {y_{i} - {A_{x} \cdot x_{i}}} \right)}}}}^{2} \cdot {w\left( {x,x_{i}} \right)}} & (5) \end{matrix}$

-   where w(x, x_(i))=e^(−(x−x) ^(i) ⁾ ^(T) ^(Λ(x−x) ^(i) ⁾ -   subject to det(Λ)=c -   c=positive constant

Salient features of the optimization problem in Equation (5) are that (i) in this new setting, Λ or really the shape matrix S as well as the regression parameter matrix A_(x) are jointly optimized; and (ii) the constraint c placed on the determinant of Λ fixes the size of the neighborhood.

Standard search-based constrained optimization techniques with a suitable choice of a starting point, can be used to determine the optimum S and A_(x).

Note finally that although the embodiment and examples have been described for linear regression, in principle the same technique readily extends for nonlinear regression. The elements of the matrix A would be simply replaced by the parameters of the nonlinear approximation function.

Revisiting the regression problem in FIG. 3, the training data was fed to both classical local linear regression as well as to the proposed regression with neighborhood shaping. FIG. 5 visualizes the results. In both cases, for the same neighborhood size (appropriately parameterized in either cost function) contours of the distance function are overlaid on the regression data. Clearly, in the proposed case, regression data points in the “linear region” afford more weight and hence the regression succeeds in perfectly approximating the function.

With reference to FIG. 1, illustrated is a flow chart of an exemplary image processing regression method incorporating neighborhood shaping as discussed hereto. By way of example, the method will be described with particular reference to a printing system, however, the image processing method is not limited to printing and can be applied to image processing in general.

In the context of a printing system, the image processing regression method illustrated in FIG. 1 generates a multidimensional color transformation associating a plurality of device dependent color space values with a plurality of respective device independent color space values.

To print a particular image, the inverse of the color transformation generated by the method of FIG. 1 is accessed to transform an image represented in device independent color space to device, i.e. printer, dependent color space for rendering/printing on the printer. As previously discussed, the device dependent color space representations of the image numerically indicate a relative amount of their corresponding printer colors, e.g. CMYK, necessary to print the original image represented in device independent color space.

To generate the printer characterization transform, computer readable instructions are executed in the following sequence:

Initially, the printer multidimensional color characterization transform generation algorithm starts 2.

Next, a set Γ of training samples (x_(i), y_(i)), 1≦i≦T, is received where x_(i) represents input color data to the multidimensional color transformation, and y_(i) represents corresponding output color data from the multidimensional color transformation 4.

Next, a parameterized form of a regression function f(x) that approximates the multidimensional color transformation is selected 6.

Next, an input color x is received 8.

Next, a cost function C representing a localized error produced by the regression function f(x) on the training set Γ is generated 10, where the localized error is a function of both the parameters of f(x) and a shaping function that defines the shape and orientation of a neighborhood of training data localized around the input color x.

Next, the parameters of the regression function f(x) and shaping function to jointly minimize the cost function C are derived 12.

Next, an output color y is generated 14 by calculating f(x) using the derived parameters of the regression function f(x) and shaping function to jointly minimize the cost function C.

Finally, the printer multidimensional color characterization transform generation algorithm ends 16.

With reference to FIG. 6, illustrated is a black diagram of an image processing system using an exemplary regression method as discussed hereto and described with reference to FIG. 6.

In operation, the printing system receives a digital input 100, represented in device independent color space, and processes 102 the device independent color space representation of the digital input image 100 to generate a pixel representation of the digital input image suitable for printing on printing device 106 to generate a hardcopy output 108 of the digital input image 100.

The image processing path 102 can reside and be executed on a DFE (Digital Front End), and/or the printing device 106. However, as will be understood by those of skill in the art, any computer related device capable of executing instructions can be used to process the image data.

As shown in FIG. 6, the image processing path includes a multidimensional color transformation, e.g. a look-up-table, which incorporates data generated by the color transformation derivation module to produce device dependent color space representations of the digital input image. Notably, the color transformation derivation module approximates a multidimensional color transformation according to the methods described in this disclosure and specifically illustrated in FIG. 5.

After the digital input image is processed by the multidimensional color transformation module to produce device dependent color space representations of the digital input image, the image data is processed according to specific tone reproduction curves 112 and halftoning algorithms 114 to generate pixel data to be rendered on the printing device 106.

Minimization of the cost function in Eqn (5) is not a trivial problem. According to another aspect of this disclosure, provided is an efficient and stable algorithmic solution to optimize the regression matrix A_(x) as well as the shape parameters in Λ such that a local minimum of the regression cost function in Eqn (5) is obtained. While simple gradient descent based optimization algorithms are available in software packages and can be applied to this problem, the merits of the disclosed algorithm, and exemplary embodiments thereof, lie in the fact that the solution exploits the structure of the regression cost function in Eqn (5).

In general, most prior optimization techniques can only hope to arrive at local minima (See C. G. Atkeson, A. W. Moore and S. Schaal, “Locally weighted learning,” Artificial Intelligence Review, 1997). One class of cost functions that is particularly convenient for optimization is the set of convex functions. A convex function defined over a convex constraint set exhibits the desirable property that a local minimum is also a global minimum within that constraint set. Furthermore, convex functions lend themselves to tractable and efficient solutions.

In solving the optimization problem, the following key observations are made. When seen as a joint function of A_(x) and Λ, the regression cost function C(A_(x), Λ) in Eqn (5) is NOT convex in its arguments. However, for a fixed choice of A_(x)=A₀, it can be shown that the function C(A₀, Λ) is a convex exponential function of Λ, and the constraint set {Λ: det(Λ)=c} is also a convex set. Likewise, for a fixed Λ=Λ₀, the function C(A_(x), Λ₀) is a convex quadratic function in terms of A_(x).

Thus, while Eqn (5) does not have joint convexity, it does have separable convexity. This observation gives rise to the following efficient optimization algorithm:

Alternating Least Squares (ALS) algorithm for optimizing regression and shape parameters:

Step 1: Initialization−Set Λ=I_(n) (the n×n identity matrix); A_(x)=A₀, any matrix in R^(m×n);

Step 2: Least Squares for A_(x): For the fixed -Λ obtained/initialized in the previous step, solve a convex least squares problem to obtain a unique solution A_(x)=A*;

Step 3: Least Squares for Λ: For optimized A_(x) in Step 2, solve a convex least squares problem to obtain a unique solution Λ=Λ*; and

Step 4: Check stopping criterion: Evaluate the gradient ∇_((p)) C(A _(x),Λ),p=vec(A _(x),Λ).

For a sufficiently small ε>0,

IF ∥∇_((p))C(A_(x),Λ)∥<ε THEN Terminate algorithm with Aopt=A*, and Λopt=Λ*

ELSE Go to Step 2.

In the aforementioned algorithm, it is implicitly assumed that the determinant of Λ is always fixed to 1—the scalar parameter Λ in Eqn (4) may be used to control the determinant.

Some salient features of the disclosed ALS algorithm include the following:

The separable convexity of the regression cost function individually is exploited to solve a series of alternating convex least squares problems—much commercial software is available to solve the individual convex optimizations problems very efficiently.

The stopping criterion of measuring the magnitude of the gradient (or its closeness to zero) is essentially the same as looking for changes in the error function value as a function of the optimization parameters but numerically more stable.

Looking for a gradient of zero implicitly assumes that local minima will be attained naturally as a result of this algorithm and forms an attractive feature of the proposed algorithm. (This fact can be proven formally.)

Standard gradient based approaches may be used to optimize p or in other words A_(x) and Λ jointly will yield solutions very sensitive to initialization, i.e. quality of minima could be poor. The disclosed approach is relatively more robust to the initial choice of Λ.

Notably, in practice, standard optimization algorithms can take a very long time to converge because they do not approach a convex problem at any stage. The disclosed algorithms exploit separable convexity to achieve significant efficiency. This is particularly advantageous in color management applications, where it is often desired to evaluate the local linear regression function for a large number of nodes of a multidimensional lookup table, where compute-intensive approaches are undesirable.

While the preferred embodiment has been described in terms of a linear regression function, the technique is applicable also to nonlinear regression which is appropriately parametrized and yields a convex optimization problem for fixed shape parameters.

With reference to FIG. 7, illustrated is a flow chart of an exemplary image processing method including the utilization of an alternating least squares algorithm as disclosed. By way of example, the method will be described with reference to a printing system, however, the image processing method is not limited to printing and can be applied to image processing in general.

In the context of a printing system, the image processing method in FIG. 7 generates a multidimensional color transformation associating a plurality of device dependent color space values with a plurality of respective device independent color space values.

To print a particular image, the color transformation generated by the method of FIG. 7 is accessed to transform an image represented in device independent color space to device, i.e. printer, dependent color space for rendering/printing on the printer. As previously discussed, the device dependent color space representations of the image numerically indicate a relative amount of their corresponding printer colors, e.g. CMYK, necessary to print the original image represented in device independent color space.

To generate the printer characterization transform, computer readable instructions are executed in the following sequence:

Initially, the printer multidimensional color characterization transform generation algorithm starts 2.

Next, at block 4, a set Γ of training samples (x_(i), y_(i)), 1≦i≦T, is received where x_(i) represents input color data to the multidimensional color transformation, and y_(i) represents corresponding output color data from the multidimensional color transformation 4.

Next, at block 8, an input color x is received.

Next, at block 200, a regression function f(x) is selected, where the regression function is parameterized by a regression matrix A_(x) that approximates the multidimensional transformation at the input data point x.

Next, at block 10, a cost function C representing a localized error produced by the regression function f(x) on the training set Γ is generated, where the localized error is a function of both the parameters of f(x) and a shaping function that defines the shape and orientation of a neighborhood of training data localized around the input color x.

Next, at block 202, the method derives the parameters of the regression function f(x) and shaping function to jointly minimize the cost function C by using an iterative alternating least squares algorithm to determine the elements of the regression matrix A_(x) and a shaping matrix Λ associated with the shaping function.

Next, at block 204, an output color y is generated by calculating f(x) using the derived parameters of the regression function f(x).

Finally, at block 16, the printer multidimensional color characterization transform generation algorithm ends.

The color characterization transform can be represented as a multidimensional color lookup table (LUT), as is commonly done in color management applications. In this case, the aforementioned process is used to evaluate the output y for each node x of the LUT lattice. Once the lookup table is constructed in this fashion, an image can be processed through the LUT using standard lookup and interpolation techniques.

It will be appreciated that various of the above-disclosed and other features and functions, or alternatives thereof, may be desirably combined into many other different systems or applications. Also that various presently unforeseen or unanticipated alternatives, modifications, variations or improvements therein may be subsequently made by those skilled in the art which are also intended to be encompassed by the following claims. 

1. A regression method for approximating a multidimensional color transformation comprising: (a) receiving a set Γ of training samples (x_(i), y_(i)), 1≦i≦T, where x_(i) represents input color data to the multidimensional color transformation, and y_(i) represents corresponding output color data from the multidimensional color transformation; (b) receiving an input color x; (c) selecting a regression function f(x) parameterized by a regression matrix A_(x) that approximates the multidimensional color transformation at the input color x; (d) generating a cost function C representing a localized error produced by the regression function f(x) on the training set Γ, where the localized error is a function of both the parameters of f(x) and a shaping function that defines the shape and orientation of a neighborhood of training data localized around the input color x; (e) deriving the parameters of the regression function f(x) and shaping function to jointly minimize the cost function C by using an iterative alternating least squares algorithm to determine the elements of the regression matrix A_(x) and the shaping function; and (f) generating an output color y by calculating f(x) using the derived parameters of step (e).
 2. The method according to claim 1, wherein the shaping function in step (d) is a shaping matrix S.
 3. The method according to claim 2, wherein the regression function in step (c) is a local linear regression function represented as y=ƒ(x)=A _(X) ·x, xεR ^(n) , yεR ^(m) , AεR ^(m×n).
 4. The method according to claim 3, wherein the cost function C in step (d) is defined as ${C = {{C\left( {A_{x},\Lambda} \right)} = {\frac{1}{T}{\sum\limits_{i = 1}^{T}{{{y_{i} - {A_{x} \cdot x_{i}}}}^{2} \cdot {w\left( {x,x_{i}} \right)}}}}}},$ where w(x, x_(i))=e^(−(x−x) ^(i) ⁾ ^(T) ^(Λ(x−x) ^(i) ⁾, Λ=λS^(T)S, S is a shaping matrix, and the minimization of C is subject to the constraint that the determinant of Λ=c, where c is a positive constant.
 5. The method of claim 1, wherein the multidimensional color transformation is a mapping from a device dependent color space to a device independent color space.
 6. The method according to claim 1, the iterative alternating least squares algorithm of step (e) comprising: (e1) setting Λ=I_(n) and A_(x)=A₀, where I_(n) is an n×n identity matrix and A₀ is any matrix in R^(m×n); (e2) for the fixed Λ of step (e1), obtaining a unique solution A_(x)=A* by utilizing a convex least squares algorithm; (e3) for the unique solution A_(x)=A* obtained in step (e2), obtaining a unique solution Λ=Λ* by utilizing a convex least squares algorithm; (e4) evaluating the gradient ∇_((p))C(A_(x),Λ),p =vec (A_(x),Λ), and IF ∥∇_((p))C(A_(x),Λ)∥≦ε THEN perform steps (e2) and (e3) again, OTHERWISE set A_(x) =A* and Λ=Λ*, where ε>0 and the determinant of Λ is fixed at
 1. 7. A non-transitory computer-readable medium storing a computer program, that when executed by a computer, causes the computer to execute a regression method for approximating a multidimensional color transformation, the method comprising: (a) receiving a set Γ of training samples (x_(i), y_(i)), 1≦i ≦T, where x_(i) represents input color data to the multidimensional color transformation, and y_(i) represents corresponding output color data from the multidimensional color transformation; (b) receiving an input color x; (c) selecting a regression function f(x) parameterized by a regression matrix A_(x) that approximates the multidimensional color transformation at the input color x; (d) generating a cost function C representing a localized error produced by the regression function f(x) on the training set Γ, where the localized error is a function of both the parameters of f(x) and a shaping function that defines the shape and orientation of a neighborhood of training data localized around the input color x; (e) deriving the parameters of the regression function f(x) and shaping function to jointly minimize the cost function C by using an iterative alternating least squares algorithm to determine the elements of the regression matrix A_(x) and the shaping function; and (f) generating an output color y by calculating f(x) using the derived parameters of step (e).
 8. The non-transitory computer-readable medium according to claim 7, wherein the shaping function in step (d) is a shaping matrix S.
 9. The non-transitory computer-readable medium according to claim 8, the regression function in step (c) is a local linear regression function represented as y=ƒ(x)=A _(X) ·x, xεR ^(n) , yεR ^(m) , AεR ^(m×n).
 10. The non-transitory computer-readable medium according to claim 9, wherein the cost function C in step (d) is defined as ${C = {{C\left( {A_{x},\Lambda} \right)} = {\frac{1}{T}{\sum\limits_{i = 1}^{T}{{{y_{i} - {A_{x} \cdot x_{i}}}}^{2} \cdot {w\left( {x,x_{i}} \right)}}}}}},$ where w(x,x_(i))=e^(−(x−x) ^(i) ⁾ ^(T) ^(Λ(x−x) ^(i) ⁾, Λ=λS^(T)S, S is a shaping matrix, and the minimization of C is subject to the constraint that the determinant of Λ=c, where c is a positive constant.
 11. The non-transitory computer-readable medium according to claim 7, wherein the multidimensional color transformation is a mapping from a device dependent color space to a device independent color space.
 12. The non-transitory computer-readable medium according to claim 7, wherein the iterative alternating least squares algorithm of step (e) comprises: (e1) setting Λ=I_(n) and A_(x) =A₀, where I_(n) is an n×n identity matrix and A₀ is any matrix in R^(m×n); (e2) for the fixed A of step (e1), obtaining a unique solution A_(x)=A* by utilizing a convex least squares algorithm; (e3) for the unique solution A_(x)=A* obtained in step (e2), obtaining a unique solution Λ=Λ* by utilizing a convex least squares algorithm; (e4) evaluating the gradient ∇_((p))C(A_(x),Λ),p=vec(A_(x), Λ), and IF ∥∇_((p))C(A_(x),Λ)∥>ε THEN perform steps (e2) and (e3) again, OTHERWISE set A_(x) =A* and Λ=Λ*, where ε>0 and the determinant of Λ is fixed at
 1. 13. An image processing method for rendering an image on an image output device comprising: receiving a device independent color space representation of the image; processing the device independent color space representation of the image through a characterization transform associated with the image output device to generate a device dependent color space representation of the image, the characterization transform associating a plurality of device independent color space values with a plurality of respective device dependent color space values, the multidimensional color transformation generated by performing a method comprising: (a) receiving a set Γ of training samples (x_(i), y_(i)), 1≦i ≦T, where x_(i) represents input color data to the multidimensional color transformation, and y_(i) represents corresponding output color data from the multidimensional color transformation; (b) receiving an input color x; (c) selecting a regression function f(x) parameterized by a regression matrix A_(x) that approximates the multidimensional color transformation at the input color x; (d) generating a cost function C representing a localized error produced by the regression function f(x) on the training set Γ, where the localized error is a function of both the parameters of f(x) and a shaping function that defines the shape and orientation of a neighborhood of training data localized around the input color x; (e) deriving the parameters of the regression function f(x) and shaping function to jointly minimize the cost function C by using an iterative alternating least squares algorithm to determine the elements of the regression matrix A_(x) and the shaping function; and (f) generating an output color y by calculating f(x) using the derived parameters of step (e).
 14. The image processing method according to claim 13, wherein the shaping function in step (d) is a shaping matrix S.
 15. The image processing method according to claim 14, wherein the regression function in step (c) is a local linear regression function represented as y=ƒ(x)=A _(X) ·x, xεR ^(n) , yεR ^(m) , AεR ^(m×n).
 16. The image processing method according to claim 13, wherein the cost function C in step (d) is defined as ${C = {{C\left( {A_{x},\Lambda} \right)} = {\frac{1}{T}{\sum\limits_{i = 1}^{T}{{{y_{i} - {A_{x} \cdot x_{i}}}}^{2} \cdot {w\left( {x,x_{i}} \right)}}}}}},$ where w(x,x_(i))=e^(−(x−x) ^(i) ⁾ ^(T) ^(Λ(x−x) ^(i) ⁾, Λ=λS^(T)S, S is a shaping matrix, and the minimization of C is subject to the constraint that the determinant of Λ=c, where c is a positive constant.
 17. The image processing method according to claim 13, wherein the multidimensional color transformation is a mapping from a device dependent color space to a device independent color space.
 18. The image processing method according to claim 13, wherein the iterative alternating least squares algorithm of step (e) comprises: (e1) setting Λ=I_(n) and A_(x)=A₀, where I_(n) is an n×n identity matrix and A₀ is any matrix in R^(m×n); (e2) for the fixed Λ of step (e1), obtaining a unique solution A_(x)=A* by utilizing a convex least squares algorithm; (e3) for the unique solution A_(x)=A* obtained in step (e2), obtaining a unique solution Λ=Λ* by utilizing a convex least squares algorithm; (e4) evaluating the gradient ∇_((p))C(A_(x),Λ), p=vec(A_(x),Λ), and IF ∥∇_((p))C(A_(x),Λ)∥>ε THEN perform steps (e2) and (e3) again, OTHERWISE set A_(x)=A* and Λ=Λ*, where ε>0 and the determinant of Λ is fixed at
 1. 19. The image processing method according to claim 13, wherein the characterization transformation is represented as a multidimensional lookup table.
 20. A non-transitory computer-readable medium storing a computer program that when executed by a computer, causes the computer to perform a color transformation for rendering an image on an image output device, the method of deriving the color transformation comprising: (a) receiving a set Γ of training samples (x_(i), y_(i)), 1≦i≦T, where x_(i) represents input color data to the multidimensional color transformation, and y_(i) represents corresponding output color data from the multidimensional color transformation; (b) receiving an input color x; (c) selecting a regression function f(x) parameterized by a regression matrix A_(x) that approximates the multidimensional color transformation at the input color x; (d) generating a cost function C representing a localized error produced by the regression function f(x) on the training set Γ, where the localized error is a function of both the parameters of f(x) and a shaping function that defines the shape and orientation of a neighborhood of training data localized around the input color x; (e) deriving the parameters of the regression function f(x) and shaping function to jointly minimize the cost function C by using an iterative alternating least squares algorithm to determine the elements of the regression matrix A_(x) and a shaping matrix Λ associated with the shaping function; and (f) generating an output color y by calculating f(x) using the derived parameters of step (e).
 21. The non-transitory computer-readable medium according to claim 20, wherein the shaping function in step (d) is a shaping matrix S.
 22. The non-transitory computer-readable medium according to claim 21, wherein the regression function in step (c) is a local linear regression function represented as y=ƒ(x)=A _(X) ·x, xεR ^(n) , yεR ^(m) , AεR ^(m×n).
 23. The non-transitory computer-readable medium according to claim 22, wherein the cost function C in step (d) is defined as ${C = {{C\left( {A_{x},\Lambda} \right)} = {\frac{1}{T}{\sum\limits_{i = 1}^{T}{{{y_{i} - {A_{x} \cdot x_{i}}}}^{2} \cdot {w\left( {x,x_{i}} \right)}}}}}},$ where w(x, x_(i))=e^(−(x−x) ^(i) ⁾ ^(T) ^(Λ(x−x) ^(i) ⁾, Λ=λS^(T)S, S is a shaping matrix, and the minimization of C is subject to the constraint that the determinant of Λ=c, where c is a positive constant.
 24. The non-transitory computer-readable medium according to claim 23, wherein the multidimensional color transformation is a mapping from a device dependent color space to a device independent color space.
 25. The non-transitory computer-readable medium according to claim 20, wherein the iterative alternating least squares algorithm of step (e) comprises: (e1) setting Λ=I_(n) and A_(x) =A₀, where I_(n) is an n×n identity matrix and A₀ is any matrix in R^(m×n); (e2) for the fixed Λ of step (e1), obtaining a unique solution A_(x) =A* by utilizing a convex least squares algorithm; (e3) for the unique solution A_(x)=A* obtained in step (e2), obtaining a unique solution Λ=Λ* by utilizing a convex least squares algorithm; (e4) evaluating the gradient ∇_((p))C(A_(x),Λ), p=vec(A_(x),Λ), and IF ∥∇_((p))C(A_(x),Λ)∥>ε THEN perform steps (e2) and (e3) again, OTHERWISE set A_(x) =A* and Λ=Λ*, where ε>0 and the determinant of Λ is fixed at
 1. 26. The non-transitory computer-readable medium according to claim 20 wherein the multidimensional color transformation is represented as a multidimensional lookup table. 